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Abstract 


We have chosen the Large Scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) code to calculate 
the coalescence of silicon carbide (SiC), silicon oxide dust (SiO) in the AGB stellar wind. LAMMPS is a classical 
molecular dynamics simulation code. At the same time, we consider the effect of temperature on the evolution of 
molecular dynamics. We also calculated the temperature change of non-spherical SiC, SiO dust coalescence. The 
condensation temperature range of SiC dust in the AGB stellar wind is [300-500]k and [900-1100]k for SiO. 
Finally, the infrared spectrum of SiC was calculated using Gaussian 16 software. The 77SiC, 70Si3C3, and 
121Si;C3 models have clear characteristic peaks of infrared spectra responding at 5, 8.6, 11.3, 15, 19, and 37 um. 
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1. Introduction 


During the asymptotic giant branch (AGB) phase, stars with 
an initial mass of up to 8 Mo lose their outer stellar layers by 
means of a stellar outflow or wind. This creates an extended 
circumstellar envelope (CSE). The CSE is a rich astrochemical 
laboratory that has detected more than 100 molecules as well as 
newly formed dust species. The chemical content of the CSE is 
determined by the elemental carbon-to-oxygen ratio (C/O) of 
the AGB star itself, with C/O < 1 leading to an oxygen-rich 
outflow and C/O > 1 leading to a carbon-rich outflow (see, 
e.g., Habing & Olofsson 2003). 

During the thermally pulsing AGB phase, low- and inter- 
mediate-mass stars develop strong mass loss and eject their 
envelope to become a short-lived (Siess et al. 2022). Stars are 
enshrouded by an optically thick dust shell, in which solids 
condense from the gas phase as tiny dust grains (Bowen 1988; 
Gail & Sedlmayr 1999). Dust is divided into O-rich dust and 
C-rich dust, pure thermal decomposition of olivine results in the 
formation of more silicon oxide dust (SiO), Mg, Fe, O (Gail & 
Sedlmayr 1999). C-rich dust mainly consists of amorphous 
carbon, and silicon carbide (SiC; Van de Sande et al. 2019). Since 
there are no condensation seeds in the gas simply outflowing from 
the star. In a series of physicochemical events, silicate dust grains 
eventually grow in the cooler outflows. The first step in the series 
of events leading to the growth of silicate dust grains in a cooling 
outflow must be the formation of seed particles. Because SiO is 
the most abundant of all gas-phase refractory elements that may 


infrared: ISM — 


(ISM:) dust — extinction 


be involved in this process, there has been much discussion about 
whether it is the clusters formed from this substance that are the 
basis of all the processes? (Gail & Sedlmayr 1986; Ali & 
Castleman 2005; Nuth & Ferguson 2006; Paquette et al. 2011; 
Goumans & Bromley 2012) 

Gail et al. (2013) show an upper limit to the stability of SiO 
based on new measurements of the vapor pressure of SiO 
molecules against solid SiO. SiO may be the first substance to 
nucleate out of a mixture of oxygen-rich elements (Nuth & 
Donn 1981). The dominant population (>0.9) of presolar silicon 
carbide (SiC) are believed to have originated from low-mass AGB 
stars with close-to-solar metallicities (Liu & Stephan 2022). Some 
90 years ago, Wildt (1933) considered the possibility that solid 
SiC grains might form in N-type stars. It is now recognized that 
SiC solids are a major dust species, second only to amorphous 
carbon particles, condensing in the cold atmospheres of mass- 
losing, carbon-rich AGB stars (Nanni et al. 2021; Chen et al. 
2022). Spectroscopic observations have revealed the widespread 
presence of SiC dusts in carbon stars, with a prominent 11.3 jum 
emission feature that characterizes the Si-C stretching of SiC 
solids (Speck et al. 1997; Mutschke et al. 1999). However, the 
exact mass fraction of condensate in the form of SiC in carbon 
stars is imprecise and depends on the stellar mass and metallicity 
(Nanni et al. 2021). Radiative transfer modeling of the observed 
infrared emission from carbon stars suggests that the mass ratio of 
SiC to amorphous carbon may be as high as 0.25 in the MW 
(Groenewegen et al. 1998), 0.43 in the LMC, and 0.11 in the 
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SMC (Groenewegen et al. 2009; Nanni et al. 2019). Theoretical 
dust yield calculations predict that the mass fraction of SiC 
produced in carbon stars out of the total dust is 0.25 of the total 
dust (i.e., silicon carbide plus amorphous carbon) for carbon stars 
with an initial mass of 3 Mọ in the vicinity of the solar metallicity 
(Nanni et al. 2013; Zhukovska & Henning 2013). Small 
molecules containing silicon and carbon are thought to be 
important components of interstellar particles. Some of them have 
been detected in the interstellar environments of late stars by 
rotational spectroscopy, such as SiC, SiC, SiC, c-SiC3, SiC4, 
while centrosymmetric species, such as C3, C4, Cs, SixCz, Si2C3, 
can only be detected by vibrational excursions (mainly in the 
infrared) (Witsch et al. 2017). Witsch et al. (2017) present 
laboratory data on the asymmetric stretching mode of SiC at 
8.6 um. They used a pulsed Nd:YAG laser to vaporize a solid Si 
target in a dilute methane sample exposed to helium buffer gas. 
Si2C was formed in the adiabatic expansion of a supersonic jet, 
and radiation from a quantum cascade laser was used to record 
rotationally resolved spectra. 

Theoretical calculations (Gail & Sedlmayr 1986) based on 
the results of Nuth & Donn (1981) show that nucleation of SiO 
should not start until the temperature drops below ~600 K, 
whereas observed condensation temperatures of periastron dust 
shells are typically in the range of 900 K and sometimes even 
higher. 

The total ice number density and number of monolayers 
throughout the outflow depend on the temperature profile of the 
dust, with clear differences between the two amorphous 
carbonaceous grains and SiC dust (Van de Sande et al. 
2019). The dust temperature profile of C-rich dust SiC is 
350 K-500 K, according to the thesis of Van de Sande et al. 
(2019) can be obtained, where the C-rich dust mainly consists 
of large amounts of amorphous carbon and SiC. The dust 
temperature profile of O-rich dust is 550-800 K. The basic case 
of primary grain condensation can be best illustrated by stellar 
wind. Here, stellar wind refers to the basic situation of stellar 
outflow. For example, dust gas temperature and material 
dilution phenomena occur. There is a monotonic decrease in 
both gas temperature and matter density, which leads to a 
recolonization of the plasma and the appearance of diatoms, 
polyatoms, triatoms, tetraatoms, pentatoms, and hexatoms as 
the distance increases. As the distance increases, diatomic 
molecules, polyatomic molecules and finally composite 
molecules are formed. At certain distances, at temperatures 
well below 1500 K, certain atoms and molecules combine into 
stable clusters. The combination of atoms and molecules into 
stable clusters can lead to the formation of homogeneous or 
heterogeneous macroscopic samples by the addition of further 
similar or dissimilar substances (SedImayr 1994). 

In recent years, many research groups have studied the dust 
and heavy elements CH, Be: He, IN, DN, and 160) yield of 
AGB stars (Ferrarotti & Gail 2006; Karakas & Lattanzio 2007; 
Ventura et al. 2020, 2021; Cinquegrana & Karakas 2022), and 
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in supernovae Wu et al. (2021). We estimated the abundance 
value of SiC in the solar system based on the spectrum, 
provided that the spectrum was corrected. Until now, only a 
few groups study the SiC and SiO dust condensation processes 
in AGB stellar winds, and infrared spectra of SiC dust. It has 
important research significance to the role of temperature on 
molecular dynamics, not only SiC, but also SiO, and the 
comparative analysis of its radiation spectrum. The study of the 
content and structure and spectra of SiC has a deeper role and 
impact. 

The pressure magnitudes are [10 *-10-°] atm and 
[10 °-10-°] atm, respectively (Sedlmayr 1994; Van de Sande 
et al. 2019; Gail et al. 2020; Siess et al. 2022). Due to the very 
low pressure, down to 1076 atm, this effect will be ignored in 
this article. The topic of SiC and SiO dust condensation 
processes is especially important and fascinating. Because we 
can better understand the origin of SiC and SiO and the basic 
laws of stars, the evolution of ISM. It can be used as the James 
Webb Space Telescope (JWST) to observe the structure and 
study of dust. In Section 2, the physical parameters of the 
model used are summarized. Detailed results are discussed in 
Section 3. The conclusions are outlined in Section 4. 


2. Models and Methods 


We choose the Large-scale Atomic/Molecular Massively 
Parallel Simulator (LAMMPS) code to compute the SiC, SiO 
dust condensation process in the AGB stellar winds. LAMMPS 
is a classical molecular dynamics (MD) simulation code 
focusing on materials modeling. It was designed to efficiently 
run on parallel computers and to be easy to extend and modify. 
For simulation, it can be viscosity, diffusion coefficient, 
energy, radial distribution function, melting and solidification 
of spherical or ellipsoidal particles (https: / /docs.lammps.org/ 
Manual.html) (Plimpton 1995; Denniston & Navid 2022). The 
model parameters are set as follows: SiC and SiO atoms were 
uniformly scattered with an initial spacing of 2 A in a 
simulation box sized 20 x 20 x 20 A? with periodic bound- 
aries. The box size of the particle heat bath is 40 x 40 x 40. 
The timestep is 0.05 fs. The velocity of the atoms conforms to a 
Gaussian distribution. The system synthesis is NPT, because 
we want to consider the condensation of SiC, SiO dust at 
different temperatures and pressures. The temperature range for 
condensation of SiC dust in the AGB stellar wind is [300—500] 
k and [900-1100]k for SiO. The numbers of atoms of SiC and 
SiO dust are 100 and 121, respectively. Regarding the potential 
function, David et al. (2012) used LAMMPS to study the 
simulation of the initial oxidation process regarding O, and 
H,0 molecules on the SiC surface. ReaxFF is able to naturally 
simulate the breaking and formation of interatomic bonds. This 
ReaxFF is also able to describe both Si-C—O and Si-O-H bond 
interactions. The ReaxFF force field, which is based on the 
principle of determining the bond levels according to the 
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distance between atoms, describes the interactions of the 
chemical system related to bonding according to the bond 
levels, and can smooth the transition of the system energies 
during bonding and bond breaking. Compared with the 
conventional force field, the ReaxFF force field provides a 
better description of the polarization because the atomic 
charges are dynamically optimized using the electronegativity 
balance method at each time step in the simulation. Compared 
with quantum mechanics, the ReaxFF molecular simulation is 
more suitable for complex systems because the reaction paths 
can be determined entirely by system driving (Du & 
Wang 2018). In most cases, LAMMPS numerically integrates 
the Newtonian equations of motion for particles such as atoms, 
molecules or macroscopic particles. At the beginning of the 
program, the expressions for the resolved force and potential 
energy in the discrete coordinate positions are computed and 
stored in an array. The simulation process uses the sample strip 
interpolation is used in the simulation to obtain accurate force 
and potential energy values. Mean square displacement (MSD): 
the square of the average distance (r) is traveled by the particles 
at time t interval. 


N 
(PO) = = a RO. (1) 
i=l 


Here, N denotes the number of atoms. Diffusion coefficient: 
(r°(t)) =2dD, where d is the dimension. Velocity autocorrela- 
tion function (VACF) (v(A)v(0)): 


D= Í =d: (2) 


Infrared radiation spectra obtained by Gaussian 16 calculation. 
A uniform broadening (Lorentzian line shape) imposed by 
Gaussian software to simulate real spectral lines. In general, 
electronic spectra (UV-Vis, ECD) are broadened with Gaussian 
functions, and vibrational spectra (IR, Raman, VCD) are 
broadened with Lorentzian functions, which decay more slowly 
than Gaussian functions. According to the uniform bandwidth 
—Lorentzian line equation: 
Tmax (=) 


(w — wo’) + (4) 


The lifetime bandwidth over the wavenumber is 
v=5.3/ t(cm~' ps_'). For a quantum state with a lifetime of 
T. When we take the frequency mean interval w — wọ? as 10. In 
Equation (3), we take half the size of radiation intensity (w), 
like [(w)/2. We can obtain Equation (4). Imax represents the 
maximum radiation intensity. About (w), the bandwidth v 
range obtained after taking the maximum and minimum values 
of the radiation intensity is 0.1—-0.5 nm. 
Tmax 1 


= 4 
7 2I (w) (w — wo?) A(w — wo) 1o 


Iw) = (3) 
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Infrared spectral (IR) calculations were performed on the final 
structure of LAMMPS after applying the temperature. IR 
spectroscopy originates from vibrational jumps and consists of 
measuring the absorption of matter exposed to light at a range 
of IR frequencies. The total amount of energy absorbed is 
called the intensity, which is plotted as the height of the peak at 
that frequency. The spectrum thus obtained is a unique 
identifier of a molecule’s constituent atoms and bonding (i.e., 
its structure). we applied the density functional theory (DFT). 
Regarding the generalization we chose B3LYP with the basis 
group 6-21g. The B3LYP function is a popular and well- 
regarded function that is primarily calculated for metal oxides 
and other inorganic compounds. The researchers calculated the 
asymmetric model also using the B3LYP function (Lee et al. 
1988). The first step calculation is to obtain the self-consistent 
energy, optimized atomic configuration and spin density. 
Required convergence on the root mean square density matrix 
was less than 1078. A single scale factor of 0.965 was used to 
compare the DFT harmonic wavelength NDFT (cm ') with the 
experimental data (Ota 2015). Corrected wavenumber N is 
obtained simply by N (cm ')=NDFT (cm') x 0.965. 
Wavelength A is obtained by À (um) = 10,000/N(em~')(Ota 
2019). 


3. Results 


Our results show that 100 atoms of SiC, SiO are stabilized 
by molecular dynamics stabilization calculations to obtain the 
final stabilized structure. Among them, the 77SiC, 70Si3;C3, and 
121Si;C3 models have characteristic peaks of infrared spectra 
responding at 5, 11.3, 15, 19, and 37 um, which are represented 
by curves with different colors respectively. 


3.1. SiC and SiO Molecular Dynamics and Spectral 
Calculations for SiC 


The classical process of dust formation is divided into two 
steps. First, critical clusters are formed that act as seeds of 
condensation, and then the dust grains are grown to macro- 
scopic size by the addition of small entities called monomers. 
Second, the dust grains are grown to macroscopic size by 
adding small monomers entities (Siess et al. 2022). In Figure 1, 
the initial structure of the 121Si;C3; models comes from the 
work of Wu et al. (2022). SizC3 was chosen because it is the 
smallest cluster with a three-dimensional structure (see 
Gobrecht et al. (2017)), so we could study the spectra of the 
low-end clusters. Our reasons for choosing the three molecular 
structures 77SiC, 70Si3C3 and 121Si3C3 as representatives for 
the IR spectral calculations include the following two main 
points. The first is that few researchers have calculated the IR 
spectra of the above three molecular structures. The second 
reason is that the above three molecules are easy to calculate 
the non-spherical spectra and easy to replicate the established 
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(a) (b) 


(c) (d) 


(e) (f) 


Figure 1. Subfigure (a) shows the initial distribution of 100SiC molecules in a box size of 40 x 40 x 40, and subfigure (b) shows the final coalescence of 77SiC 
molecules into nuclei in the box. Subfigure (c) shows the initial distribution of 121Si3;C3 molecules in a box size of 40 x 40 x 40, and subfigure (d) shows the final 
coalescence of 121Si3C3 molecules into multiple small dust nuclei in the box. Here the number 121 indicates the total number of atoms of the structure. Similarly, 
subfigure (e) shows the initial distribution of 100SiO molecules in a box size of 40 x 40 x 40, and subfigure (f) shows the final coalescence of 100SiO molecules into 
multiple small dust nuclei in the box. In the figure, the yellow atoms are C, the red atoms are O, and the blue atoms are Si. 
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Figure 2. The left and right graphs show the trend of temperature with time for SiC, SiO molecules. The black line represents SiC molecules, and the red line 


represents SiO molecules. 


model, which is convenient for us to calculate other physical 
properties at a later stage. So we chose the above three 
molecular models as the main model calculation. For 
LAMMPS, the basic principle of molecular dynamics is the 
use of Newtonian mechanical equations to describe the 
positions, velocities, and interatomic interactions of atoms, 
and by solving these equations, it is possible to simulate the 
kinetic behavior of molecules and to compute the kinetic 
behavior of large-size systems. We can see that the atoms 
coalesce into small SiC and SiO clusters due to interactions. It 
mainly consists of the attractive and repulsive forces between 
molecules. In the gaseous state, the smaller the intermolecular 
distance, the weaker the intermolecular interaction, and 
vice versa, the stronger the interaction force. When the 
molecules are at high temperature or high pressure, the 
intermolecular spacing gradually shrinks, which will produce 
a large amount of absorbed or released heat, the gravitational 
force is formed by the charge interaction of molecules, when 
the distance between molecules is not changing, and finally 
leads to changes in morphology, the formation of a solid. In 
Figure 2, the temperature is calculated by dividing the kinetic 
energy by a certain number of degrees of freedom (and the 
Boltzmann constant). Since kinetic energy is a function of the 
particle velocity, it is often necessary to distinguish between the 
advective velocity of a particle (due to the aggregation motion 
of some particles) and its thermal velocity. The sum of the two 
is the total velocity of the particle, but the latter is usually the 
temperature that we want to calculate. Let us perform a 
constant NPT system to update the position, velocity, and 
angular velocity of a finite-size x-spherical particle for each 
time step (refer tohttps://docs.lammps.org). The initial 


temperature of the SiC molecule drops from 500 to 300K, 
where the pressure is so small that we ignore its effect (refer to 
the work of Hanine et al. (2020), and the highest temperature in 
the figure can even be as high as 800K, due to the 
intermolecular forces being so strong that they cause the 
temperature to rise. Similarly, the initial temperature of the SiO 
molecule decreases from 1100 to 900K, with a maximum 
temperature of 1600 K. In Figure 3, the total energy of the two 
molecules is demonstrated to reach a stabilization value, 
showing the stabilization of the final structure and the 
molecular dynamics reaching an equilibrium state. Both models 
are down to equilibrium trend. For example, the final total 
energy of the SiC dust cluster is as low as 
—1.2 x 10*kcal mole !, the SiC cluster mass size 
3.31 x 107” kg. SiO cluster mass size ranges from 
3.65 x 107” to 7.3 x 10°** kg. Energy includes kinetic energy 
E1=0.5mv’, potential energy E2=mgh, where “m” is the 
mass of the atom and “v” is the velocity value. “h” is the 
distance traveled by the atom. Rotational energy E3 = 0.51”. 
Here, “E” represents the rotational kinetic energy of the atom, 
“P is the rotational moment of inertia of the atom, “w” is the 
angular velocity of the atom. 

In Figure 4, the horizontal coordinate indicates the 
wavelength and the vertical coordinate indicates the magnitude 
of the radiation intensity. 77SiC model has peaks at 
wavelengths of 6, 8.6, 10, 15, 20, 31, 33 um, respectively. In 
particular, there is a strong characteristic peak at 15 um, which 
we find to be due to a large amount of Si—C bond stretching. As 
shown in Figure 5 right, the Si-C, Si-Si bond length scaling of 
the structure on the right side of the blue solid line is more 
significant than that on the left side of the scale. 121Si3C3 
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Figure 3. The evolution curve of the total energy (E) of SiC, and SiO molecules 
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Figure 4. The infrared spectra of model 77SiC are represented by the red curve, and the remaining models 121Si;C3 and 70Si3C3 are represented by green and dark 
blue, respectively. The figure in the lower right corner shows a comparison of the infrared spectra of the three models. 


model has peaks at wavelengths of 8.6, 12, 16, 19, 23, 27, 29, 
37, and 66 um, respectively. In particular, there are two strong 
characteristic peaks at 19, 37 um, which we consider to be due 
to the large amount of Si-C, Si-Si bond stretching. As shown 


in Figure 5 left, the Si-C, C=C bond length scaling of the 
structure on the right side of the blue solid line is more 
significant than that on the left side of the scale. The 70Si3C3 
model has peaks at wavelengths of 5, 11.3, 14, 17, 18, 21, 29, 
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Figure 5. The picture on the left shows model 121Si;C; and the picture on the right shows model 77SiC, where the gray atoms are represented as C atoms and the dark 


green atoms are Si atoms. 


53, 60 um, respectively. In particular, there are two strong 
characteristic peaks at 5, 11.3 wm, which we consider to be due 
to a large amount of Si-C bond stretching. Among the three 
models, the 70Si;C3 model has the highest radiation intensity. 
121Si3C3 and 77SiC models have decreasing radiation intensity 
in that order. Regarding the radiative signature of SiC at 
11.3 um (Speck et al. 1997; Mutschke et al. 1999), only model 
70Si3zC3 of our model has a significantly radiative signature, 
and the remaining two models are apparently not significant, 
and we consider that this is due to the weaker Si-C bond 
scaling and the stronger Si-Si, C=C bond scaling. Witsch et al. 
(2017) presented laboratory results for the asymmetric stretch- 
ing mode of SiC at 8.6 um. Here we obtain this peaks for both 
7T7SiC and 121Si3;C3 models. This result provided theoretical 
support for later observations of dust in the ISM. We consider 
still the peak due to Si-C bond stretching. This calculation is of 
great significance for the future budget of the output of stellar 
dust, galaxy dust, and the influence of dust radiation spectra on 
observations. For example, it can serve as a reference for dust 
observations by the JWST, the Large Sky Area Multi—Object 
Fiber Spectroscopy Telescope (LAMOST), and the upcoming 
China Space Station Optical Survey Telescope (CSST). 

We estimated the abundance value of SiC in the solar system 
based on the spectrum and calculated it as follows: the infrared 
emissivity per unit mass of SiC molecules is 


a/m = EDan). © 
The mass absorption coefficient is 
Cabs 3 Cabs 
Tap z aap’ . 
The absorption cross-section is 
Cavs = Ha?’ Qavs (A). (7) 


The above formula can obtain: 


G /m) = 22A 
ap 


IIB, (T). (8) 
If SiC has a radiation flux of Fy at 11.3 um, the distance 
between a star and Earth is d. 


_ Msic(ijy/™) 


E 
Ald? 


(9) 
According to formulas (8) and (9), the following equation can 
be obtained: 


4ap 
3Qans(A)B\(T) 


Here, we have made corrections to F, by multiplying it by 
10-°° to match the observed spectrum. Using the distance d 
between the Sun and the Earth (d= 1.47 x 10''m). The 
frequency is 2.75 x 10'* Hz. The absorption efficiency factor 
Qavs is 0.7, dust density p is 3.22 gem >, dust radius a is 
1.17 x 107° m, and the Planck function B,(T) is 2.28 x 107". 
By substituting the above parameters into Equation (10), the 
abundance value of SiC (Msic) can be obtained. Msic is 
4.75 x 1075, which happens to be within the corresponding 
range of 10 °-10-? for Ferrarotti & Gail (2006) job. 


Msic = Kd? (10) 


4. Conclusions 


Through MD simulations, we modeled the process of atoms 
coalescing into small clusters of SiC and SiO in AGB stellar 
winds. It was observed that the initial temperature of SiC 
molecules decreased from 500 to 300K, while that of SiO 
molecules decreased from 1100 to 900K. SiC molecules 
condensed entirely into a single cluster, whereas SiO molecules 
condensed into multiple small clusters within an evolution time 
of up to 10° fs. 
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Furthermore, we analyzed the infrared spectrum of SiC 
nanoparticles. It was found that the 77SiC, 70Si;C3, and 
121Siz;C3 molecules exhibited characteristic peaks at 5, 8.6, 
11.3, 15, 19, and 37 um. Specifically, the 77SiC molecule 
exhibited peaks at 8.6 and 15 um, the 70Si3;C3 molecule 
exhibited peaks at 5 and 11.3 wm, and the 121Siz;C; molecule 
exhibited peaks at 8.6, 19, and 37 wm. Many of these peaks 
corresponded to the Si-C bond stretching mode. 
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